function [t, phi, Vt, It] = memristor(a, b, A, F, phi0, tspan)
% myMemristor 使用自定义的 RK4 方法计算忆阻器在正弦输入下的状态演化
%
% [t, phi, Vt, It] = myMemristor(a, b, A, F, phi0, tspan)
%
% 输入参数:
%   a, b   : 忆阻器导纳函数 W(φ)= a + b*tanh(φ) 的参数
%   A      : 正弦信号幅值
%   F      : 正弦信号频率 (Hz)
%   phi0   : 状态变量 φ 的初始值
%   tspan  : 求解时间区间, 例如 [0, 2*(1/F)]
%
% 输出参数:
%   t      : 时间向量
%   phi    : 状态变量 φ(t)
%   Vt     : 输入电压 V(t)=A*sin(2*pi*F*t)
%   It     : 输出电流 I(t) = (a + b*tanh(φ))*V(t)

    % 定义输入电压函数
    V_fun = @(tt) A .* sin(2*pi*F*tt);

    % 定义状态方程 dphi/dt，使用 RK4 求解
    % 状态方程： dphi/dt = -0.2*phi^3 + V(t) + phi
    odefun = @(tt, phiVar) -0.2*phiVar.^3 + V_fun(tt) + phiVar;

    % 使用自定义 RK4 求解器计算 φ(t)
    [t, phi] = rk4(odefun, tspan, phi0);

    % 计算对应时刻的输入电压 V(t) 和输出电流 I(t)
    Vt = V_fun(t);
    It = (a + b*tanh(phi)).*Vt;
end
